clear all
set more off
cd /Users/zimaoxiao/Dropbox/CC_Yield_Predict/replication_package/  // Navigate to replication folder on your own machine

cap mkdir output

qui {
* prepare for mapping
use county_fips yield2020 ccimpact_corn_585_lin trends_corn_585_lin ccimpact_corn_585_quad trends_corn_585_quad using data/figure/baseline_2050, clear
merge 1:1 county_fips using data/geo/cty_eastern_data, nogenerate
foreach i in yield2020 ccimpact_corn_585_lin trends_corn_585_lin ccimpact_corn_585_quad trends_corn_585_quad {
	replace `i' = -9999 if `i' == .
}

* map: 2020 yields (in bu/acre)
grmap yield2020 using data/geo/US_County_LowRes_2013coord ///
, id(_ID) osize(vvthin ...) ocolor("189 189 189" ...) ///
	clmethod(custom) clbreaks(-9999 50 70 90 110 130 150 170 190 210 230) /// 
	fcolor(white ///
	"255 255 229" "247 252 185" "217 240 163" "173 221 142" "120 198 121" "65 171 93" "35 132 67" "0 104 55" "0 69 41") ///
	leg(order(11 "210 -- 230" 10 "190 -- 210" 9 "170 -- 190" 8 "150 -- 170" 7 "130 -- 150" ///
	6 "110 -- 130" 5 "90 -- 110" 4 "70 -- 90" 3 "50 -- 70" 2 "Not in sample")) ///			
	line(data("data/geo/st_eastern_coord") size(vthin ...)) ///
	graphr(margin(zero)) ///
	legend(size(*1.2) position(3) ring(1)) ///
	title("(a) 2020 corn yields in bu/acre", size(*0.8))
graph save "Graph" "output/FigureS4(a).gph", replace
	
* map: CC effects (linear trends)
grmap ccimpact_corn_585_lin using "data/geo/US_County_LowRes_2013coord" ///
	, id(_ID) osize(vvthin ...) ocolor("189 189 189" ...) ///
	clmethod(custom) clbreaks(-9999 -1.2 -1 -0.8 -0.6 -0.4 -0.2 -0.1 -0.05 0 0.05 0.1) /// 
	fcolor(white ///
	"103 0 13" "165 15 21" "203 24 29" "239 59 44" "251 106 74" "252 146 114" "252 187 161" "254 224 210" ///
	"247 252 245" "229 245 224") ///
	leg(order(12 "+5% -- +10%" 11 "0% -- +5%" 10 "-5% -- 0%" ///
	9 "-10% -- -5%" 8 "-20% -- -10%" 7 "-40% -- -20%" 6 "-60% -- -40%" ///
	5 "-80% -- -60%" 4 "-100% -- -80%" 3 "-120% -- -100%" ///
	2 "Not in sample")) ///			
	line(data("data/geo/st_eastern_coord") size(vthin ...)) ///
	graphr(margin(zero)) ///
	legend(size(*1.2) position(3) ring(1)) ///
	title("(b) CC effects (linear trends)", size(*0.8))
graph save "Graph" "output/FigureS4(b).gph", replace

* map: CC + Tech effects (linear trends)
grmap trends_corn_585_lin using "data/geo/US_County_LowRes_2013coord" ///
	, id(_ID) osize(vvthin ...) ocolor("189 189 189" ...) ///
	clmethod(custom) clbreaks(-9999 -0.5 -0.3 -0.1 -0.05 0 0.05 0.1 0.3 0.5 0.7 0.9) /// 
	fcolor(white ///
	"252 146 114" "252 187 161" "254 224 210" "255 245 240" ///
	"247 252 245" "229 245 224" "199 233 192" "161 217 155" "116 196 118" "65 171 93") ///
	leg(order(12 "+70% -- +90%" 11 "+50% -- +70%" 10 "+30% -- +50%" 9 "+10% -- +30%" ///
	8 "+5% -- +10%" 7 "0% -- +5%" 6 "-5% -- 0%" 5 "-10% -- -5%" 4 "-30% -- -10%" 3 "-50% -- -30%" ///
	2 "Not in sample")) ///			
	line(data("data/geo/st_eastern_coord") size(vthin ...)) ///
	graphr(margin(zero)) ///
	legend(size(*1.2) position(3) ring(1)) ///
	title("(c) CC + Tech effects (linear trends)", size(*0.8))
graph save "Graph" "output/FigureS4(c).gph", replace

* map: CC effects (quadratic trends)
grmap ccimpact_corn_585_quad using "data/geo/US_County_LowRes_2013coord" ///
	, id(_ID) osize(vvthin ...) ocolor("189 189 189" ...) ///
	clmethod(custom) clbreaks(-9999 -1.2 -1 -0.8 -0.6 -0.4 -0.2 -0.1 -0.05 0 0.05 0.1) /// 
	fcolor(white ///
	"165 15 21" "203 24 29" "239 59 44" "251 106 74" "252 146 114" "252 187 161" "254 224 210" "255 245 240" ///
	"247 252 245" "229 245 224") ///
	leg(order(12 "+5% -- +10%" 11 "0% -- +5%" 10 "-5% -- 0%" 9 "-10% -- -5%" ///
	8 "-20% -- -10%" 7 "-40% -- -20%" 6 "-60% -- -40%" 5 "-80% -- -60%" 4 "-100% -- -80%" 3 "-120% -- -100%" ///
	2 "Not in sample")) ///			
	line(data("data/geo/st_eastern_coord") size(vthin ...)) ///
	graphr(margin(zero)) ///
	legend(size(*1.2) position(3) ring(1)) ///
	title("(d) CC effects (quadratic trends)", size(*0.8))
graph save "Graph" "output/FigureS4(d).gph", replace
	
* map: CC + Tech effects (linear trends)
grmap trends_corn_585_quad using "data/geo/US_County_LowRes_2013coord" ///
	, id(_ID) osize(vvthin ...) ocolor("189 189 189" ...) ///
	clmethod(custom) clbreaks(-9999 -1.5 -1.3 -1.1 -0.9 -0.7 -0.5 -0.3 -0.1 -0.05 0 0.05 0.1 0.3 0.5 0.7 0.9) /// 
	fcolor(white ///
	"103 0 13" "165 15 21" "203 24 29" "239 59 44" "251 106 74" "252 146 114" "252 187 161" "254 224 210" "255 245 240" ///
	"247 252 245" "229 245 224" "199 233 192" "161 217 155" "116 196 118" "65 171 93") ///
	leg(order(17 "+70% -- +90%" 16 "+50% -- +70%" 15 "+30% -- +50%" 14 "+10% -- +30%" 13 "+5% -- +10%" 12 "0% -- +5%" ///
	11 "-5% -- 0%" 10 "-10% -- -5%" 9 "-30% -- -10%" 8 "-50% -- -30%" 7 "-70% -- -50%" 6 "-90% -- -70%" ///
	5 "-110% -- -90%" 4 "-130% -- -110%" 3 "-150% -- -130%" ///
	2 "Not in sample")) ///			
	line(data("data/geo/st_eastern_coord") size(vthin ...)) ///
	graphr(margin(zero)) ///
	legend(size(*1.2) position(3) ring(1)) ///
	title("(e) CC + Tech effects (linear trends)", size(*0.8))
graph save "Graph" "output/FigureS4(e).gph", replace		

* ------------------------------
* combine to get Figure S4
gr combine ///
output/FigureS4(a).gph ///
output/FigureS4(b).gph ///
output/FigureS4(c).gph ///
output/FigureS4(d).gph ///
output/FigureS4(e).gph ///
, row(2) col(3)
noi graph export "output/FigureS4.png", as(png) name("Graph")

erase "output/FigureS4(a).gph"
erase "output/FigureS4(b).gph"
erase "output/FigureS4(c).gph"
erase "output/FigureS4(d).gph"
erase "output/FigureS4(e).gph"
}
	
*** EOF
